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The crystallization conjecture: a review 


Xavier Blanc and Mathieu Lewin 


In this article we describe the crystallization conjecture. It states that, in appropriate physical conditions, 
interacting particles always place themselves into periodic configurations, breaking thereby the natural 
translation-invariance of the system. This famous problem is still largely open. Mathematically, it amounts 
to studying the minima of a real-valued function defined on R 3JV where N is the number of particles, which 
tends to infinity. We review the existing literature and mention several related open problems, of which 
many have not been thoroughly studied. 
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At the microscopic scale, most crystals are composed of atoms which are arranged on a periodic 
lattice. This specific geometric structure has important consequences at the macroscopic scale. 
For instance, in snowflakes the atoms are arranged on an hexagonal lattice, which explains the 
beautiful six-pointed figures that can be found in nature. The aim of crystallography is to study 
those periodic structures and their properties at larger scales. 

2014 was declared the year of crystallography by UNESCO [236] and this gives us the opportunity 
to draw attention to a difficult mathematical conjecture, also important from a physical point of 
view, which has been studied a lot without being completely solved. While crystallographers study 
the properties of some periodic arrangements and compare them, there remains a more fundamental 
question: why is it favorable (at low temperature) for the atoms to spontaneously arrange themselves 
on a periodic array? This periodic order seems to only appear in the limit of a large number of 
particles, which makes the question particularly difficult. 

In this article we rigorously formulate this long-standing problem and we make a review of the 
existing results as well as of the remaining open questions. We will mostly discuss the simplest model 
(classical particles interacting with a two-body potential at zero temperature), before addressing 
more advanced situations (for instance quantum systems and/or positive temperature). 


1. The classical model 

1.1. Energy. Let us consider a set of N classical identical particles in (in practice d = 1,2,3), 
interacting by pairs through a potential V depending only on the distance between them. We 
denote by x\,..., Xn G K d and pi,... ,pn G the positions and momenta of these particles. The 
model to be used is that of the Hamiltonian dynamics, based on the energy 

N I 12 

U N (xi, . . . ,X N ,P!,. . .,p N ) = ^2 + X! V(\Xi-Xj I). (1) 

2=1 1 <i<j<N 

Here m is the mass of the particles and | • | is the Euclidean norm of R d . 

At zero temperature, the equilibrium states are the minima of "Hat, which all satisfy pi = ■ ■ ■ = 
Pn = 0. If one is only interested in those, it is therefore sufficient to consider the potential energy 

£n(x X N )= ^2 V(\Xi-Xj\), 

l<i<j<N 

and to understand how the x^s solving the minimization problem 

E(N) = infj£jv(a;i, • ■ - ,xn), xi ,... ,x N G R d j, (2) 

are arranged in W l in the limit N —> +oo. Let us note that En is invariant under translations and 
rotations. Any configuration may be rotated and translated by a fixed vector without changing the 
total energy. Minimizers of © are thus not unique. At positive temperature one should consider 
the Gibbs measure exp (—T-Ln/T), as will be discussed later in Section [TT~T1 below. 

In practice the potential V depends on the type of atoms and is not explicitly known. As atoms 
are not elementary particles, V cannot be deduced from first principles. It is therefore important 
to obtain mathematical results which are sufficiently generic with regards to V. 
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Figure 1. Left: the Lennard-Jones potential ©• Right: a minimizer for the variational problem m, 
computed numerically in [20] , with N = 100 and d = 2. The particles seem to arrange themselves on an 
hexagonal lattice, and to form a large cluster having the shape of an hexagon. 

Qualitatively, the function V is usually assumed to be positive (repulsive) at small distances 
and negative (attractive) at large distances. Since the interaction between two atoms which are 
far from each other is small, we assume that V(r) —> 0 as r —> +oo. A typical and very popular 
example is the Lennard-Jones potential 




(3) 


drawn in Figure [I] The behavior at infinity in r -6 mimics the Van der Waals interaction, that 
is, the one for radially symmetric neutral particles. The behavior at r = 0 is, on the other hand, 
completely empirical. The number ro > 0 is the equilibrium distance for two isolated particles. 
It may be seen from Figure [T] that for this specific potential in the plane, the solutions Xt to the 
minimization problem (©) are approximately located on an hexagonal lattice and that they moreover 
form a big cluster having the shape of an hexagon. Nobody knows how to prove these observations 
rigorously and this is essentially the crystallization conjecture that will be discussed in this article. 

For the rest of the article, we consider a general radial potential V that tends to zero at infinity. 
Some assumptions are however necessary to ensure that our question is well posed. 

1.2. Binding: existence of minimizers for E(N). The first assumption that we need is 
that particles can bind, that is, E(N) should possess at least one minimizer (maybe only for a 
subsequence Nj —> oo). For instance, if V >0, then E(N) = e = 0 for any N £ N, but the 
minimization problem @ has no solution. Indeed, the inhmum of £n(xi, ..., Xn) is reached only 
when the distances between particles Xi tend to infinity. It is therefore mandatory to assume that 
min V < 0 (or to force the particles to stay together, as will be discussed in Section [Till . In physical 
situations the interaction is attractive at some distance and we therefore always have 

min V < 0. 

Note that E( 2) = min V, hence this implies that E(2) < 0. 

Next we remark that the energy is subadditive, that is, it satisfies 


\/N,P> 1 , E(N + P) < E(N) + E(P). 


( 4 ) 
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This inequality is shown by sending P particles at infinity and using the fact that V —> 0. In other 
words, we write 


E(N + P) < lim £ n+p (xi, ...,x N ,yi + r,y P + r) = £ N (x i, ...,x N ) + £p{yi ,..., yp) 

| r | —>-oo 


and get 0 after optimizing with respect to the Xj’s and ijj's. Applying 0 inductively, we find 
that 

E(N) < [N/2\E{2) <(N- 1)^^. 

The optimal energy E(N) is negative and bounded above by a term which behaves linearly in the 
particle number N. Since we require that min V < 0, we deduce that 


E(N) < E(2) 

N ~ 2 



< 0 . 


(5) 


When the strict binding inequalities 


E(N) < E(K) + E(N - K) 


( 6 ) 


hold for all K = 1,..., TV — 1 and when V is a continuous function on (0,oo), then E(N) can be 
shown to possess at least one minimizerQ The verification of ([6]) can be complicated for a general 
V , but it is easy if V < 0 at infinity. In this case, © can be proved by induction on N , using 
that if two groups of particles are far away, they always attract each other. We conclude from this 
discussion that when V is continuous on (0, oo) and negative at infinity, then E(N) always possesses 
minimizers and satisfies (0. 

It can also be proved that is non-decreasing, which gives us an inequality in the reverse 

order: 


E(N) 

N 


>(N- 1 ) 


E( 2) 

2 


(7) 


The latter is indeed obvious from the formula of £/v, since each of the N(N — l)/2 terms in the 
sum can be bounded from below by V{\ Xi — Xj\) > min(V) = E(2). As we will explain in the next 
section, the lower bound dTJ) is not optimal in physically interesting cases. 


1.3. Stability and the behavior of E(N) for large N. We have seen that E(N) is bounded 
above by a linear term in N. This linear behavior is indeed the interesting physical case and we 
will always require that the following limit 


^OO 


lim 

N— )-+oo N 


( 8 ) 


exists and is finite. The reason is the following: if we gather two macroscopic identical systems (for 
a “real life” object, N « 10 23 ), the formation energy is by definition equal to to 2E(N) — E(2N) > 0, 
which may be arbitrarily large if E(N ) is super-linear Q Note that eoo < 0 from 0, since we require 
that min(V) < 0. 

x The idea is to consider a minimizing sequence and to study whether some particles escape. If all the Xj s are 
uniformly bounded (after applying an appropriate translation) then we get the existence of a minimizer using the 
continuity of V. If K particles escape, then we get a contradiction from GQ. 

2 For instance 2 E(N) - E(2N) ~ C(2“ - 2 )N a if E(N) ~ -CN a with a > 1 and C > 0. 
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It is clear that the existence of the limit © implies that there is a lower bound of the form 


E(N) > -CN, 


( 9 ) 


but the converse assertion is actually true and it is sometimes called Fekete’s subadditive lemma [85]. 
Using that E(N) < |_ N/K\E(K ) for every fixed K and taking the limit N —> oo, we even deduce 


^-oo 


inf 

N>1 


E(N) 
N 


In other words, C = —e x is the optimal constant in ©. So we have to restrict ourselves to the 
potentials V for which © is satisfied, and those are called stable in the literature. The lower 
bound © can be rewritten in the form 


£ V(\xi-Xj\) > -CN (10) 

for all N and all xi,...,xn £ Stable potentials have been widely studied since the 60s ©Hi 

[Ml HEl ECO H© H© HZ] . 

The simplest example of a stable potential is 


V = Vi + V 2 , with Vi > 0, V 2 > 0 and 


V 2 < + 00 , 


/R d 


( 11 ) 


where V 2 denotes the Fourier transform of x >->■ V 2 (|x|), cf. [209] Prop. 3.2.7]. For V 2 , the proof 
relies on the observation that 


N N 


v^Xi-Xji) = ~ y F (°) 


l<i<j<N 


i=l j =1 

2(27r) d / 2 J Rd 


N 


V 2 (k) \Ye 
3 =1 


2 n N 

lk ' Xi dk- — U(0) > -yU(0). (12) 


However, there are many physical potentials that cannot be written in the form dill) . 

Another famous example is that of a potential which is non-integrable at 0 but is bounded from 
below by an integrable function at infinity. Namely, V is stable when 


-C 
~V2 (r) 


for 0 < r < a 
for a < r < b 
for r > b 


with ip 1 and f 2 two positive decreasing functions on (0, a) and (b, 00 ) respectively, such that 

pa poo 

/ f\{r) r d_1 dr = 00 , / f 2 (r)r d ~ 1 dr<oo, 

Jo Jb 

see [65] and [591 App. A]. The idea behind these conditions is that the fast increase at zero prevents 
the particles from being too close to each other in average. The integrability of V at infinity 
then makes the double sum in £n behave like N. The Lennard-Jones potential © satisfies these 
conditions in dimensions d < 5, and it is therefore stable. 
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1.4. Formation of a macroscopic object. The assumptions © © above imply the existence 
of a thermodynamic limit ([8]), but they do not ensure that a macroscopic object is formed in this 
limit. Their aim is actually to avoid a collapse of the system by preventing particles to be too close 
to each other. It is still possible that, in the optimal configuration, particles do not stay close. This 
situation should not be allowed. 

In order to have the formation of a macroscopic object, we want that the minimizing configura¬ 
tion of N particles fill a volume of size N, in the limit N —> + 00 . Moreover, the particles should be 
“evenly spaced” in this volume, as it is clear in the example of the Lennard-Jones potential shown 
in Figure [TJ 

The mathematical formulation of this property is not unique. One possibility is to apply a 
dilation of factor N~ 1 ^ d to the optimal configuration (this is a way to pass to the macroscopic 
scale), and to consider the empirical measure 



*=1 


(13) 


We ask if, after extracting a subsequence Nj and translating the whole system by a vector t 3 , 




M weakly-* in the sense of measures, with M £ L°°(K d ), 


supp(M) compact and 


M(x) dx 


1. (14) 


This means that the macroscopic object is included in the support of M, and that, at this scale the 
system is continuous, with the function M as local density. If (fTdl) is satisfied, a macroscopic object 
has formed. The knowledge about the positions of particles is very crude, because of the dilation 
of factor N~ 1/d , which does not account for the local behavior of the system. 

Finding conditions on the potential which imply the existence of the weak limit m is an 
important problem. However, it has never been, to our knowledge, studied mathematically. For 
the Lennard-Jones potential Vlj in dimension 2, Figure [l] indicates that the limit M is proportional 
to the characteristic function of a hexagon. In this example, the shape of the support of M, which 
is visible at the macroscopic scale (as for instance snow flake structure), is a manifestation of the 
crystalline order at the microscopic scale. 

In the next section we discuss the crystallization conjecture, which concerns the microscopic 
properties of the system. We will come back to the macroscopic scale later in Section 12.81 below. 


2. The crystallization conjecture 

2.1. Formulation. We now come to the question which has been intensively studied since the 
1960s, without being solved [ 2351 T96j . This question is: does the system become periodic in the 
limit N —> 00 ? This may be formulated as follows. Let us denote by 

N 

hn = ^2 ( 15 ) 

i=l 

the empirical measure associated with the solution xi ,..., Xn of problem ({2]) . Note that, contrary 
to (11311 . we do not use any dilation, and this means that we study the system at the microscopic 
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scale. We ask if, after extracting a subsequence Nj and translating the whole system by a vector 
T b 

VNji'-Tj) ^ H (16) 

locally, where /i is a locally finite measure. We say that crystallization occurs if y is periodic, that 
is, if there exists a discrete subgroup 


G = 


d 

, j=i 


Tin G Z 


c R d , 


(17) 


generated by d independent vectors v±, ...,Vd & K. d , such that p(- + g) = H f° r all 9 G G. In order 
to avoid trivial cases, we assume here that G is the maximal group satisfying this property. Put 
differently, the period is supposed to be minimal. The invariance under the action of G does not 
imply that the particles are located on the vertices of a periodic lattice. This would correspond to 
the stronger hypothesis 

A 4 = 'y 1 fig+m (18) 

g eG 

for some fixed vector y £ R d , defining the position of the lattice in space. For instance, it is possible 
to have 3 particles in the unit cell of the lattice G, which are repeated periodically, as in Figure [2] 
In such a case, the configuration of particles is the superposition of 3 shifted crystalline lattices, 
and the measure has the form 


^ = £< 

g&G 


g+y 


-£< 

seG 


'g+y+T + £1 bg+y+r’- 
g&G 




Figure 2. Example of a periodic configuration in 2D. 

In the special case where the particles are exactly on the nodes of the lattice G, as in (fT51) . we 
use the word Bravais lattice or mono-atomic lattice. For instance, in dimension 3, the simple cubic 
lattice (SC), face-centered cubic lattice (FCC) and body-centered cubic lattice (BCC) are all Bravais 
lattices. On the other hand, the hexagonal close packed lattice (HCP) is not. It is the superposition 
of two shifted Bravais lattices (Figure [3]). This configuration is the one used to pile up oranges in 
markets. 

Crystallization may be seen as a symmetry breaking of the system: the invariance of the system 
under affine isometries is lost in the process. If the positions of the particles form a periodic lattice, 
then applying a translation, rotation or reflexion to the system does not change its energy. Hence, 
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hexagonal close packed (HCP) 

Figure 3. Most common configurations in 3D. © Wikipedia, GNU FDL, by C. Dang Ngoc Chan. 


the set of minimizing lattices has the structure of the compact group (M d /G) x 0<j(M). Choosing a 
special minimizer for the positions x, at finite iV, it is possible to select one of the limiting lattices. 

The microscopic scale convergence (1161) does not give any information, in principle, about the 
behavior at the macroscopic scale, such as the convergence of the dilated measure My defined 
by d- Conversely, the convergence of M/v does not give any clue about that of //jy However, 
one actually expects that the two phenomena are related. The understanding of the link between 
these two scales is still incomplete, as we will discuss in Section 12.81 below. 

Should crystallization be proved, the next question is to know which periodic configurations 
are present in the limit (that is, what is the group G ). Another question is to know if /j has the 
particular form (1181) corresponding to a Bravais lattice. If not, one would ask how many particles 
are present in each periodic cell, and what are their positions. In physical systems, lattices with 
larger symmetry groups seem to be more common (142| . These lattices are the hexagonal and 
square lattices in 2D, and the lattices presented in Figure [3] in 3D. 

The ubiquity of crystals (at low temperature) indicates that crystallization is a universal phe¬ 
nomenon, which should occur for a wide class of interaction potentials V. As we will see, several 
mathematical works prove crystallization, but they are based on restrictive assumptions on V. To 
date, no generic class of potentials has been identified, for which crystallization can be proved. 

2.2. Convergence of /xjv and the minimal distance between the particles. Before dis¬ 
cussing existing results about the crystallization problem in itself, we discuss when the sequence 
/jjv can be shown to have subsequences that converge locally, weakly-* in the sense of measures, 
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as we required in (1161) . We need to show that the number of particles in an arbitrary fixed ball 
(centered at tn and of radius R ) is uniformly bounded in the limit N — > oo: 

#{|Xj — t n | < R} < C{R) < oo. (19) 

If this property is satisfied, by usual compactness arguments we can then construct a subsequence 
of hn{- ~ tn) which converges locally to some /i. We expect that C(R ) will behave like the volume 
\Br\ of the ball for large R, but this is not needed at this point. 

The condition dl9l) immediately follows with C(R) = |Sij|/|f3 e | + 0(R d ~ 1 ) if we can prove that 
the smallest distance between the particles does not tend to zero in the limit IV —> oo: 

min I Xi — xA > 2e > 0 

1 

for a minimizer of £n- The idea is often that a very fast blow up of V at zero should prevent that 
the particles get too close, but proving this very natural property can actually be a difficult task. 
It has been the object of several works: the Lennard-Jones potential was covered in [25111171 :1 2501 
125] 121911253 ] and more general potentials were addressed in [ 16811241112421 [3j. As we will see in 
Section EH the situation is much easier at positive temperature where it is sufficient to estimate 
the probability that some particles get very close. 

2.3. Crystallization results and sphere packing. 

In dimension d — 1. In dimension one, the problem of crystallization is rather well understood. 
The first results are due to Ventevogel and Nijboer [2381123911240| : they prove that the limit eoo 
is reached by equidistant configurations. This property is proved for a wide class of potentials 
V (they are assumed to be non-increasing up to a distance tq > 0, and non-decreasing for r > 
ro, with additional hypotheses on V"), which includes the Lennard-Jones potential Vlj. The 
convergence (1161) is not proved in these works and is still an open problem. It has been proved in 
the special case of Vlj by Gardner and Radin | |l05 j. 

For some explicit examples of potentials V (non-increasing up to a distance ro > 0, and non¬ 
decreasing for r > ro), it has been proved that the optimal configuration does not converge to a 
Bravais lattice. The limit can be clusters of particles which are globally periodic [238]. With an 
oscillating potential V, it is even possible to obtain configurations which have no periodicity [124| . 
In the latter, it is also proved that such aperiodic configurations can be found as minimizers of a 
potential V which is an arbitrarily small perturbation of a potential for which crystallization occurs. 
This indicates that crystallization is an unstable property if V is perturbed by a small but highly 
oscillating function. Thus, the conditions on V ensuring crystallization are probably complex and 
have not been completely understood yet, even in one dimension. It is commonly assumed that the 
interaction potential is smooth, stable, non-increasing up to a distance r o, and non-decreasing for 
r > ro- However, no crystallization result has been proved under these assumptions only, even in 
one dimension. 

In dimension d > 2. In higher dimensions, the problem is far from being understood. Most 
results are based on geometrical arguments, which allow to reduce the question to the sphere 
packing problem. This question consists in finding the position of non-overlapping spheres of equal 
radii giving the largest possible density. In two dimensions, the solution is precisely the hexagonal 
lattice (see Figure [4]). Thue has given two proofs of this result (in 1892 and in 1910), which both 
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happened to contain flaws. A correct proof was then provided in 1940 by Toth [2001 1BTJ . In 
dimension three, the problem is significantly more difficult. Kepler formulated it in 1611, and it is 
therefore often called Kepler’s conjecture. A computer-assisted proof was given by Hales in 1999, 
then published in 2005 in [123]. Only recently (August 2014) has it been fully validated, after eleven 
years of work by the Fly speck team [91], who managed to give a formal proof based on the softwares 
Isabell and HOL Light. An important difference with the two-dimensional case is that, in 3D, the 
problem has many solutions, including the hexagonal close packed lattice, the face centered cubic 
lattice and even non-periodic arrangements. The fact that FCC is the unique minimizer among 
Bravais lattices was proved by Gauss [ 106] . 



Figure 4. Packing of identical disks, maximizing the density. The centers of the disks lie on a hexagonal 
lattice. 


The link between the crystallization problem and the sphere packing problem has been high¬ 
lighted by Heitmann and Radin in [129 ], Indeed, if the interaction potential V is given by 

{ -Too if 0 < r < 1, 

-1 if r = 1, (20) 

0 if r > 1, 

then the particles can be considered as hard spheres of radius 1/2. These spheres tend to touch due 
to the condition F(l) = —1. The crystallization problem is thus equivalent to the sphere packing, 
and one obtains that the solution is the hexagonal lattice in 2D, and either FCC or the other sphere 
packing solutions in 3D. 

Subsequent works aimed at generalizing this result to potentials which are similar to (1201) . but 
are closer to physically realistic interactions. For instance, in [194] . Radin considered a potential 
satisfying (® for r G [0,1], which is non-decreasing for r > 1, and tends to 0 fast enough as 
r —> Too. In a famous article [229j . Theil dealt with smoother, more realistic potentials (which 
look like Vlj), in dimension two. However, he still used restrictive hypotheses on V. This work has 
been extended to dimension three recently in [90] . in which an additional three-body term is added, 
which favors particular angles between bonds. A similar strategy has been used in dimension d = 2 
in [751 fTMl [TT0] . where the optimal lattice may be a square lattice. One can therefore consider 
that the problem is not completely understood in dimension two, and completely open in dimension 
three. 

All these results in dimensions two and three rely heavily on the similarity with the sphere 
packing problem. However, it is not clear if this should be the correct physical explanation. This 
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would exclude, for instance, configurations which are periodic but not mono-atomic. In such a case, 
particles form small groups which are repeated periodically. The crystallization conjecture for a 
more general class of potentials is still an open problem. 

The sphere packing problem becomes more complex as the dimension increases. Note however 
that, in dimensions d = 4, 8, 24, there are special lattices which are believed to solve the best 
packing problem [57]. Although the sphere packing problem in high dimension plays an important 
role in information theory, it is natural to restrict ourselves to the (physically relevant) cases of 
dimension d = 1,2, 3. Indeed it has been conjectured that crystallization only occurs in small space 
dimensions for physical interactions [2241 [2341 !233| . 

2.4. A variant: minimization at fixed density. It is possible to consider a potential V which 
does not allow for the formation of a macroscopic object, if we force the particles to stay together. 
The idea is to minimize the energy while keeping the density of particles p fixed. This may be done 
by confining the particles in a large domain f l and imposing that their number be TV ~ p|f2|, where 
| Q | is the volume of Q. This way, we get a family of problems depending on the parameter p. 

To be more precise, we consider the minimization problem for TV particles in the domain 0 

Eq(N) = mf{s N (xi,... ,x N ), xi,...,x N € fij, (21) 


and we study the limit 


e(p) 


lim 

N—> oo 
\n N |— >-oo 
N/\Qn |—>p 


EnAN) 

TV 


( 22 ) 


where p > 0 is fixed and is a sequence of domains which covers the whole space in the limit 
TV —> oo. This limit should not depend on the chosen sequence. For instance, it is often assumed 
that the measure of the boundary of fljy is a lower order term compared to its volume \£In\ [209] . 
To fix the ideas, one can think of Sljv as a cube of side length ( N/p ) 1 / d , or a convex symmetric 
domain of unit volume, dilated by a factor (TV/p) 1 / 1 *, as for instance a ball of radius proportional 
to (TV/p) 1 /'*. 

Since Eq(N) > E(N) > eoo TV, it is clear that e(p) > eoo for all p > 0, where eoo is the constant 
defined in ©. In order for the limit (1^1) to exist, we also need an upper bound. Since the particles 
are confined to the domain f2, we cannot send them to infinity anymore, hence we loose most of the 
properties of Section fT72l The energy is not necessarily subadditive as in and we may well have 
E(N) > 0. Furthermore, if the function V tends to zero too slowly at infinity, each of the particles 
in the domain U interacts with many of the other particles and the energy might grow faster than 
TV. A natural condition is to assume that V is integrable at infinity: 



dr < oo. 


We can then easily find a position of the particles that will give an energy of order A, hence 
E(N) < CN. Under this assumption and the usual stability condition ®, the limit in (E^l) can be 
shown to exist and to be independent of the sequence fljsf- Potentials that are not integrable are 
also sometimes considered, but then one should divide the energy by the appropriate power of TV. 

Whether the energy behaves linearly or not, the problem is to study the behavior of the particle 
positions X\ , ...,xn solution to the minimization problem (1211) . and the questions are similar to the 
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preceding case. A difference is that the model is no more invariant under affine isometries. Differ¬ 
ent extraction of the sequence of minimizers may in principle give limiting lattices with different 
positions. In practice, the position of the limiting lattice is often determined by the choice of a 
particular sequence Djv (see [145] for a discussion of this aspect in dimension d = 1). 

Since a cube of volume N/p can be placed in a cube £l' N of volume N/p' when p' < p, 
an optimal position of the particles for f2jv can be used in f V N , and we conclude that p H > e(p) is 
non-decreasing. Assume now that V has an attractive part and that the problem E(N) studied 
before in the whole space lead to a crystal. This crystal has an average density p = n/\Q\ where 
Q is the unit cell and n is the number of particles in it (n = 1 for a Bravais lattice). Then, if we 
impose a density p < p the domain will typically be larger than the natural size of the system and 
the particles will solve the exact problem E(N), that is, we have Eq(N) = E(N). We conclude 
that e(p) = e(p) = eoo is constant for 0 < p < p. The energy only starts to grow for p > p. 

The fact that one can consider a repulsive potential changes the physical meaning of the problem. 
In particular, the relation with the sphere packing problem is less clear. In the case of the preceding 
section, it is natural to consider that the particles are attracted to each other, and behave like hard 
spheres at short range, therefore trying to maximize the density of the system. Doing so, they tend 
to maximize the number of their neighbors. Here, particles can repel each other fiercely, and tend 
to maximize their mutual distance, while staying in the domain il. Experiments and numerical 
simulations indicate, however, that here again crystallization occurs. 

In dimension one, Ventevogel and Nijboer have proved crystallization for any density p > 0 in 
the case of non-negative non-increasing convex potential [238 ]. In [23911240] . they prove the same 
result for the potential V(x) = exp(— ax 2 ) and V(x) = (/? + x 2 ) -1 , for a, 1 3 > 0, still in dimension 
one. This allows to generalize the result to any convex combination of these potential, such as 
V(x) = / 0 + °° e-“ x dp(a), for any non-negative measure p. Such potentials may be non-convex. 
In addition, they give a necessary condition for crystallization, in any dimension: if crystallization 
occurs for sufficiently small densities, and V is continuous, then V > 0. As before, the situation is 
much less clear in dimension d > 2. 


2.5. Optimal lattices and special functions. If crystallization is assumed , it is possible to 
determine the most favorable periodic configurations by comparing their energy per particle eoo 
and e(p), defined by ([8j) and ([22]) . In some cases, this question may be related to a problem in 
analytic number theory, involving special functions. 

Indeed, if the particles lie on the vertices of a Bravais lattice G (a discrete subgroup of R d such 
as (ED, with unit cell Q), the limit energy per particle reads: 

\ E V ^)- (23) 

" seG\{o} 

Finding the optimal configuration amounts to minimize this expression with respect to G. There 
is no additional constraint on G for eoo- In contrast, one needs to fix the volume of the unit cell Q 
of G to \Q\ = 1/p when the density is fixed, as for instance in the case of a repulsive potential. 


Epstein zeta function. With a Lennard-Jones type potential 


V{r) 


1 

a 



1 (ro_\ b 
b\r ) ’ 


(24) 
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where a > b > d, we get 


2 

<?eG\{o} 


Cd(g,o) Cd(g,fc) 

a b 


(25) 


where S' is a symmetric positive definite matrix of size d, which is related to the Gram matrix of 
the basis ( w j)i<j<,/> an d such that r^S 1 ^ 2 !^ = G. Here, 


Cd{S,s) 


- y 

zez d \{o} 


1 

( z T Sz ) s / 2 


(26) 


is the Epstein zeta function [8l]. Still assuming that we have crystallization on a Bravais lattice, 
the minimal energy for Vlj reads 


e LJ,C 


min 

S—S T >o 


( US, a) 


Cd(S,b) 


(27) 


Except in dimension d = 1, the solution to this problem is still unknown, even for the physically 
relevant cases a = 12 and 6 = 6. 

If the density p > 0 is fixed as discussed in Section 12.41 it is possible to consider a repulsive 
potential V(r) = r~ s with s > d. Hence, we need to minimize the value of the zeta function (1261) . 
with respect to G (that is, with respect to the matrix S ) 


e c (p,s)= min S,s). 
S=S T >o 
det(S)=p -2 


(28) 


Here, det(S) is the volume of the unit cell of the lattice to the power 2. Applying a dilation of the 
lattice, one easily proves that 

Up, s ) = p s/d U M) 

and that it is sufficient to study the problem in which the unit cell has a volume equal to 1. Without 
loss of generality, we can thus assume that det(S) = 1. There is a link with the sphere packing 
since, in the limit s —> oo, the optimal lattice converges to a solution to the d-dimensional sphere 
packing problem 21(1 . It should be noted that Cd(S. s) is not bounded. If the smallest eigenvalue 
of S reaches 0, then Q(S, s ) tends to + 00 . 

The function s <—> Cd(S, s) has an analytic continuation to the set C \ {d}. This extension has 
a simple pole at s = d, with a residue equal to 7T d / 2 r(d/2)” 1 (if det(S') = 1), and satisfies the 
functional equation 

C d (^, S )=^- d/2 ^^C d (^- 1 ,d-s) (29) 

where S' -1 is the matrix associated with the lattice G* = [k £ : k ■ g £ Z, Wg £ G} , called the 

dual lattice of G 133 ESI- Thus, it is also possible to study the minimization problem (1281) even if 
0 < s < d. As we will see below, this problem is of great importance from a physical point of view, 
particularly if d = 3 and s = 1. Formula (1291) implies that if S' is a solution to the minimization 
problem e^(l,s), then S -1 is a solution to e^(l,d— s). 

Going back to the case of the Lennard-Jones potential m, we see that, after dilating the 
problem with fixed density p > 0, it amounts to minimize the function 


eLj,c(p) = P a/d 


min 

det(S)=l 


( US, a) i^ US,b) \ 

[a P b )■ 


(30) 
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Since b > a is assumed, for large p the problem reduces to the minimization of (//(.S', a). For 
the minimization problem (1271) with no density constraint, one finally needs to consider p which 
minimizes the function p > eLj x(p)- I n this case, the optimal lattice is unknown. 

One can derive a representation of (/d as a series with exponentially decaying coefficients. The 
most widely used method is due to Ewald [83] i31] 115( and relies on the integral representation 


1 

rv*S 


1 

2r(s/2) 



r*/ 2 " 1 dr. 


For s > d, we have, if det(S) = 1, 


(31) 


US,s) 


where 


W 2 f 1 1 

r(s/2) \ s — d s 


\ J ((e d (S, t) -1 y / 2 - 1 + {e d {s-\r) - l)^- 1 ) drj , 

(32) 


O d (S,a)= J2 e~™ zTSz 

zE 1 d 


(33) 


is the Jacobi theta function. Here, we have used Poisson’s summation formula 


9d{S, a) 


a d /2 9d 



(34) 


Formula (l32l) is also meaningful for 0 < s < d and can be used to prove that Q has an analytic 
extension to C \ {d} (T has a pole at the origin which compensates for the divergent term 1/s), as 
we already mentioned. Formula (1321) is widely used by physicists. It allows to compute numerically 
the values of , s) very accurately, allowing to formulate conjectures on what should be proved. 

We are now going to describe what is expected for the minimization of the Epstein zeta function. 


Results for and 9 in dimension 2. In dimension d = 2, it has been proved by Rankin (198), 
Cassels [48], Ennola d and Diananda [65], that the hexagonal lattice is the unique minimizer of 
the zeta function, for any s > 0, when the density is fixed. In other words, we have 

C2(S,s)-C 2 (Shex,s)>0 (35) 


for all s > 0 and all S such that det(S') = 1, where 


'5'liex 


A ( 1 1 A 

W 2 1 ) 


corresponds to the hexagonal lattice. In addition, the inequality (1551) is strict if S ^ Shex, U P to 
the invariances of the problem (rotation and change of basis of the lattice). Another proof is given 
in m\- Inequality (l35l) is still valid for s = 2 where both functions have a simple pole with equal 
residue. When s —> 0, we have a divergence which needs to be dealt with, but the result is still 
true EH- We have made some numerical calculations that confirm these results, see Figure [5] 
below. 

A famous result due to Montgomery [180 1 deals with the case of a Gaussian repulsive interaction 
in dimension d = 2. In this case, the problem reduces to the study of the Jacobi theta function (1331) . 
As for the zeta function, Montgomery proves that 82(8, a) > 02(<S'hex, a) for all a > 0 and all S such 
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that det(S') = 1 (cf. Figure O. If the potential V is a positive linear combination of Gaussians, 
Montgomery’s result implies that the optimal lattice is the hexagonal one, for a fixed unit cell 
volume. For instance, using the integral formula (1321) . one recovers the previously mentioned result 
on zeta functions. 

Subtracting two Epstein zeta functions gives a function that can be expressed as an integral of 
the function O 2 (S,a) multiplied by a weight. This weight is non-negative when p is large enough. 
Using this argument, Betermin and Zhang [22] have proved that, at high density, the optimum 
of (1501) is reached by the hexagonal lattice in 2D, for the Lennard-Jones potential Vlj. Imposing 
that p is large means that particles are close to each other. Therefore, their interaction is dominated 
by the repulsive part r ~ 12 of Vlj- The energy is close to ( 2 ( 8 , 12), which, as a function of S, reaches 
its minimum for the hexagonal lattice only. On the contrary, they prove that, when p —> 0, the 
hexagonal lattice cannot be the global minimizer. For instance, if 

3 C2U2, 6 ) - C2(<S'hex,6) 

P C2(/2,12)-C2(Shex,12)’ 

the square lattice (S CBI = I 2 ) has an energy which is smaller than that of the hexagonal lattice. 
If no symmetry breaking occurs, then the square lattice becomes the minimizer. Recall, however, 
that the true e(p) is expected to be constant when p < ~p (the density of the optimal lattice at 
which (1301) is minimal). It will not coincide with eLj x(p) i J1 ^101) at small p. For a more recent work 
in the same spirit, see Ell- 

Results and conjectures for £ and Q in dimension d > 3. In dimension d > 3, some authors 
have studied the critical points and the (local or global) minima of the Jacobi theta function and 
the Epstein zeta function. In a famous article im Sarnak and Strombergsson determined special 
local minima in dimensions 4, 8 and 24 (see also [62) (58;, [63]). In dimension 3, Ennola has proved 
that the face centered cubic (FCC) lattice is a non-degenerate local minimum of £3 (5, s) for all 
s > 0 [80]. Formula (1201) implies that its dual, the BCC lattice, is also a non-degenerate local 
minimum for 0 < s < 3. In addition, based on the sphere packing problem obtained in the limit 
s —> 00, it has been shown in [2l0j that FCC is the unique global minimizer for s large enough. As 
opposed to what Ennola conjectured in [80], FCC cannot be the unique minimizer for all s > 0. 
Indeed, formula (1291) would imply that its dual, BCC, is a minimizer for some values of s. Hence, 
a more likely conjecture would be that FCC is the unique minimizer for s > 3/2, whereas BCC 
is for 0 < s < 3/2 [2181 section 5]0 If it is assumed that the minimizer has a high-symmetry 
group, and if we only compare the energies of SC, FCC and BCC, this conjecture is corroborated 
by numerical computations presented in Figure 0 It is a very important conjecture: its proof would 
be an important advance both in analytic number theory and in solid-state physics. One of the 
difficulties in the proof is that the values of the zeta function for BCC and FCC are very close to 
each other. This implies that a quantitative argument needs to be very precise. 

Similar questions may be asked about theta function (1331) . but it seems that the corresponding 
literature is far less important. In 3D, the conjecture is that FCC is the unique minimizer for any 
a > 1, whereas BCC is for a < 1 [2181 section 5]. Here again, this conjecture is confirmed by 
numerical simulations presented in Figure [5] Note that, contrary to dimension 2, the conjecture 
for the theta function does not seem to imply it for the zeta function: formula (1311) always involves 
both the lattice and its dual for different values of a. 

3 Note that BCC cannot be a non-degenerate local minimum for all s > 0 186111191 . hence it has to be a degenerate 
critical point of C/ for some s > 3. 
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Most works consider only mono-atomic lattices. This excludes the HCP (Hexagonal Close 
Packed) lattice in dimension 3, since it is not a Bravais lattice. We refer to [ 167] for an explicit 
link between zeta functions and quantum field theory, to [226| for the link with optimal quadrature 
point repartition, and to [186 . 21711218] for the link with the optimization of the determinant of 
the Laplace operator: det(—A) = e - ^ 5,0 -*. 

As a conclusion, determining the optimal periodic lattice can, for some simple potentials, be 
related to the study of special functions, the conjecture is that the minimizer can be either the 
FCC lattice, or the BCC one. This is still an open problem (in most cases), even though research 
is very active on this subject. 


2.6. Coulomb potential and Wigner crystallization. The Coulomb potential (the real inter¬ 
action between charged particles) is not covered by any of the previous results. In dimension d , the 
Coulomb potential Vj^coui by definition the Green function of the Laplace operator, that is, the 
solution to 

-AP d ,Coul = l^-Vo, 

in the sense of distribution, where |<S' d_1 | is the volume of the sphere in dimension d, 5 ' d_1 = {x £ 
R d : |a:| = 1}. We thus have 

^L.Coui(^) = -|z|, V 2 ,Coui(a;) = ^ log |x|, V dt Gou\(x) = , E for d > 3. 

\x\ 

The function V^coui is never integrable at infinity and it is not stable in dimension d < 2. In 
dimension d > 3, V^coui is non-negative, hence no macroscopic object can be formed. If we 
perform the thermodynamic limit as in Section 1271 then E(Qn,N) behaves as jV 1+2 / d . 

Physically, a macroscopic system of charged particles is never seen in vacuum. Two alternatives 
are usually considered. The first is to put them in a trap, that is, to add an external potential 
Y^jLi Kxt ( x j ) to the energy, with I4xt(a0 —> +oo when |cc| —> oo. This is now done in the labora¬ 
tory [551 11281 12501 1231] . although a large number of particles is still difficult to reach. We discuss this 
possibility in Section 12.7.11 below. Another realistic macroscopic system of electrons is when they 
are placed in an external potential describing a background of opposite charge, which compensates 
the charges of the electrons and allows for an equilibrium to form. In a metal, the background is 
composed of the positively-charged nuclei (and possibly another set of electrons that do not move). 
These other particles can be fixed or optimized (see Section 15721) . In [ 246 ]. Wigner considered the 
simpler situation of a uniform background of density p > 0, and this is the object of this section. 
This is the so-called Jellium model, in which the background is a kind of “jelly” slowing down the 
movements of the particles. As we will explain in Section 12.7.11 below, trapped systems can be 
shown to locally behave like Jellium hence, in spite of its apparent simplicity, the Jellium problem 
is of high physical relevance. 

For a general potential H, the Wigner minimization problem reads: 


EnA N ) 


= inf 
X1,...,Xn£Q 


E 

l<i<j<N 


N r 

v (\xi - Xj\) - V (\xi — y\) dy 

A _ 1 J 



V(\x - y\) dx dy >, (36) 
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with the same limit as before 

ejeii(p) = lim ( 37 ) 

N—too iV 

|Q/v |—>oo 
N/\n N \^ P 

The second term of the energy in (1551) accounts for the interaction of our N particles with the 
homogeneous background. This is a new term compared to the preceding cases. The last term is 
the energy of the background, which is constant with respect to the positions of the particles. We 
keep it in order to have a finite limit (33 (in the case we deal with here, f n V{\ x ~ y \) dxdy 
grows like N 1+2 / d ). In the limit (1571) one imposes that lV/|flvr| —> p, which means that the particle 
density is equal to that of the background. This allows to reach an electrostatic equilibrium between 
particle repulsion and the attraction of the background. One could in principle minimize over the 
domain f 1 while fixing |ff| = Np, and then use this domain Qn, but it is often assumed that Sl^v 
has a shape which is fixed (a cube or a ball for instance), and is then dilated. The limit (1571) should 
not depend on the chosen sequence Qn. 

For a stable integrable potential V, adding the background should only affect the points close to 
the boundary of and we expect that in the thermodynamic limit ej e ii(p) = e(p)—p/2 f Rd V. The 
situation is different for potentials that decay slowly. As mentioned above, Wigner was originally 
interested in the electrostatic interaction between electrons, that is, the Coulomb potential V = 
Vd,Coui- Nevertheless, the problem (1551) makes sense as soon as V behaves like |x| -s with s > d — 2 
at infinity. Screening is the main effect that will make this possible. The idea is that each particle 
will feel an effective potential that decays faster than |a;| _s , due to the cancellations induced by 
the two additional terms. The potential is expected to decay like \x\~ s ~ 2 or even |a;| -s ~ 3 if the 
configuration of the particles has sufficiently many symmetries. 

In [246], Wigner has conjectured crystallization for V = Vd.Coui, at least if 1 < d < 3 and p 
is small enough. He also suggested that in 3D, the electrons form a body centered cubic lattice 
(BCC). In 2D, the particles are expected to form a hexagonal lattice. The same is expected for 
V(|a;|) = |x| _s with d — 2 < s < d. Numerical simulations and formal computations corroborate 
Wigner’s conjecture. However, a rigorous proof is still missing in dimensions d = 2,3. 

Dimension d = 1 is simpler and has been solved by Kunz in 1974 for small densities jl4~5 j. This 
result has been generalized to any density by Aizenman and Martin [?]. At temperature T = 0, the 
particles form a lattice of step 1/p, as in the preceding sections. If T > 0, it has also been proved 
in [7] that the particle density is periodic of period 1/p. A different proof, which applies to the 
quantum case, has been proposed by Brascamp and Lieb in [39] . It is an application of their study 
of the optimality of Gaussians in some functional inequalities. 

As in Section [2751 if crystallization is assumed (on a Bravais lattice G), it is possible to compute 
the corresponding energy per unit volume. In the present case, we have 

ejeii ,g{p) = \ X! ~ [ V{x)dx + ^ [ [ V(x-y) dxdy (38) 

2 seGUQ} Jq zJqJq 

where W is the twice-screened potential 

w{x) = v(x) - 2 ] f Q v (x-v) d v + f Q f Q v ( x + y ~ dydz 

with Q the unit cell of the lattice G, satisfying \Q\ = 1/p. As expected, the two terms have the 
same effect as a Taylor expansion of V at infinity, which increases its decay and can make the series 
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in (1551) convergent (depending on the symmetry properties of the unit cell Q). If V(x) = |x| -s at 
infinity with s > d — 2, then W{x) ~ |a:| _s_2 and no assumption on Q is necessary. For s = d — 2 
(Coulomb potential), it is sufficient to choose for Q a set which is symmetric with respect to the 
origin 0, which is always possible. Doing so, W behaves like \x\^ d ~ 1 at infinity. If Q has sufficient 
symmetry properties, it is possible to express the energy with the simple screened potential 

W(x) = V (x) - |^i V {x - y) dy. 

We refer to [ 1541 App. B] for the details. It is always useful to choose for Q the Wigner-Seitz cell, 
which has the same symmetries as the lattice G m- 

When the above series converges and V(x) = |x| s , one can prove that the energy (1551) is equal 
to the analytic extension of the first term in m, that is, 


ejeii.G(p) = C d(S, s), 


(39) 


for d > 2 and d — 2 < s < d. The proof is based on the same arguments as [1541 App. B], and on 
results from [351EU GS]. The point s = d—2 which is the one interesting for applications is always on 
the left of the pole s = d, where the Zeta function is defined through analytic continuation, showing 
the importance of studying it in this region. For d = 3, the numerical simulations presented in 
this article indicate that the body centered cubic lattice (BCC) is the unique global minimizer, as 
conjectured by Wigner. Proving this fact is still an open problem. In dimension d = 2, the energy 
has a logarithmic singularity which needs to be removed, but the problem is similar. In [2l4 j. 
Sandier and Serfaty used Montgomery’s result to prove that the optimal lattice in 2D must be the 
hexagonal lattice, in the limit s —> 0. 

Physicists usually rely on the integral representation (1511) to compute the value of the zeta 
function and compare different lattices. As an example, the energies are approximately equal to 


{ 1.41865... for the simple cubic lattice (SC), 

1.44415... for the face centered cubic lattice (FCC), 

1.44423... for the body centered cubic lattice (BCC), 

in dimension d = 3 mm- Using an argument due to Onsager [155] , Lieb and Narnhofer managed 
to prove in [ 160 ] that the true energy defined by (1571l satisfies ejeii(p) > — p 1 / 3 1.4508... for any p > 0. 
This value is very close to the expected one. However, the proof of Wigner crystallization is still 
an open problem, in dimension d > 2. 

Let us point out that Wigner’s model has been recently studied and reformulated in 216 . 2031 
ESDI- In these articles, the energy ej e ii(p) is called renormalized energy and is defined directly 
on sets of infinitely many points (which need not be on a periodic lattice), without using the 
thermodynamic limit N —> oo. 


2.7. Occurrence of the crystallization problem in other situations. In this section, we 
present a few questions that reduce to the crystallization problem stated above, or to Wigner’s 
problem. This shows that these questions are rather universal. 


2.7.1. Confined systems in the mean-field limit. The crystallization problem appears in 
some dense trapped systems that are now realized in the laboratory. Here, a change of scale is 
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needed to recover a problem set in the whole space. The prototypical situation is to minimize the 
energy 

^Vext {N) = min | v< y\ x i~ x j\) + Y Ve ^ Xi ) | > ( 40 ) 

X ' \ V 1 <i<j<N i=1 I 

where I4xt is a confining potential, which tends to +oo at infinity. The coefficient 1/N multiplying 
the interaction allows for both terms to be of the same order of magnitude in the limit N —> +oo. 
This is called the mean field regime. An example of interaction V is given by V(x ) = |x| _s with 
0 < s < d, or V ( x ) = — log |x|. Most commonly used confining potentials are the harmonic potential 
V ext (x) = \x\ 2 , and the potential 


I'ext (^) 


0 if x G M , 
Too otherwise, 


which amounts to impose that all the particles stay in the bounded set A4 C The set M. can 
be a bounded domain like a ball, or a zero-measure set such as a sub-manifold of R d , of dimension 
strictly smaller than d. 

If d = 3 and A4 = S 2 is the unit sphere, and if V(x ) = |a:| _1 , this is called Thomson’s 
problem E32I- Finding the optimal positions of the particles on the sphere, even for a fixed value 
of TV, is a famous problem which has been solved only for some values of N. Many numerical 
studies of the problem have been made, giving some insight on what the optimal configurations 
should look like. This problem is related to one of the eighteen open problems mentioned by 
Smale in 1998 [2251 [TUI 118] . It naturally occurs in many different situations: it is related to the 
construction of a set of points which discretizes the sphere as uniformly as possible (the so-called 
elliptic Fekete points |211| ): in biology, this problem can explain the form of some viruses, and the 
repartition of pores on pollen grains; it is also studied in link with “colloidosomes” |66j (Figure [9]). 
If V{x) = |x| _s , with d > s, and if Af is a sub-manifold without boundary, of dimension d— 1, the 
problem is usually called Riesz problem. We refer to (126 ] for a general presentation of the problem 
and numerical simulations. Coulomb crystals in traps are now produced in the laboratory, since 

the 90s [511121IMH23H- 

For the model (140L a second-order expansion is needed to find the crystallization problem. 
Indeed, the leading order is, under appropriate assumptions, given by the mean-field theory 

limM = ^ pr mf iiity {IJj V(\x — y\) da(x) da(y) + J V ext {x) dcr(x)j := a. (41) 

measure on |A 

The measure a, a solution to this variational problem, is in general absolutely continuous with 
respect to the Lebesgue measure, and gives the average repartition of the points in space. More 
specifically, 

1 N 

< 42 > 

i= 1 

weakly-* in the sense of measures. Since in R d the points Xi will have a typical distance of order 
_/V _1 / d , the measure a plays the same role as the macroscopic measure M discussed in (fl3l) (Tl4l) . 
For instance, in the case of the Thomson problem for which V ext confines the particles to the unit 
sphere S 2 and V(|x|) = |x| -1 , the solution is unique and equal to the uniform measure a = (47r) _1 
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on S 2 . This means that the particles tend to be uniformly distributed on the sphere. For a general 
set M or general confining potential I4 xt , proofs of (14T1) and (1421) are given in [55l 1146111791 l45l l46l 
1T401PHT! 1202) . 

A change of scale is needed to be able to study accurately how the particles are organized 
at the microscopic scale. In order to do so, it is better if the potential V behaves appropriately 
under dilations. In general, one assumes that V(|x|) = |x| _s (or V(|a;|) = — log |rc| , which formally 
corresponds to the case s = 0). When I4xt is smooth, after a dilation of N~ 1 / d around a given 
point x £ R d , the problem happens to coincide with Wigner’s model in dimension d, with the local 
density p = a(x). The total limit energy is the superposition of these local problems, and one finds 

E v ext (N) = aN + ejeii(l) [ a(x) 1+ ^ dx + o(N%), (43) 

Js. d 

where a is the constant given by (ED, and where ej e ii(l) is the Jellium energy (1371) for p = 1 with 
interaction V(x) = |a:| _s (this expansion is modified in the case V{x) = — log|x|). This result has 
been recently proved by |2l6l 12X51 l201 j in the case V{x) = — log |x| in dimensions d = 1, 2, in ( 203) 
for the Coulomb potential s = d — 2, and in [l lOOj for d — 2 < s < d. 

In the case of a sub-manifold M C l d , the scaling is modified, and instead one applies a 
dilation of 7V _1 / d where d! is the dimension of A4. It is expected that the same kind of results 
hold ( 2431124411116 . 1971 [4011001 XTj , although it has not been proved yet, except in the case of the 
sphere in dimension d = 2 with V{x) = — log|x| |44) . 

The asymptotics O(N) in (1411) is only valid if V is locally integrable, so that the right side is 
finite. Several authors have studied the case of a potential which is not locally integrable, typically 
V(x) = for s > d. In the case of a submanifold A4 of dimension d !, the corresponding energy 
behaves like N s / d (or N log N for s = d'). If s > d', it was proved in [140) 1120) il26 , 147 . 1721102) 
that the corresponding term reads 

= (44 > 

where e(l) is now the minimal energy (l22l) for the problem on the whole space with V(x) = |ar|: 


e(l) = lim 
N—>oo 
|—>oo 

As we pointed out in Section 12751 the conjecture is that the particles are located on a hexagonal 
lattice in dimension d! = 2 and FCC when d! = 3. In such a case, the right side is equal to Cd'{S, s) 
with S corresponding to the optimal lattice [40j. 

Except in dimension 1 (21511148] for which the problem is better understood, it seems that none 
of these works provide any new information on the crystallization conjecture itself. Nevertheless, 
they give an important insight on how it naturally arises in many different situations. 



Xi — Xj | s , Xi £ fl 


N 


2.7.2. Vortices and crystallization in dimension 2. In dimension d = 2, the crystallization 
problem appears when studying fast rotating Bose-Einstein condensates or superconductors in 
large magnetic field. Vortices are created, and their number grows with the rotation speed (or the 
magnetic field intensity). When this number becomes large, they seem to form a hexagonal lattice, 
called Abrikosov lattice in this context [2]. 
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In fast rotating Bose-Einstein condensates, vortices may be modeled as classical particles in¬ 
teracting via a potential. The corresponding energy may be computed using the Jacobi theta 
function (l33l) [ 6 ]. In this context, Montgomery’s result explains why the vortices should form a 
hexagonal lattice (see Figure [Toll. 

Vortex patterns for the Ginzburg-Landau equation of superconductivity have been widely stud¬ 
ied in the mathematics literature (see [23] [25, the first articles on the subject, using simplified 
models). Under some constraints on the magnetic field, there is a finite number of vortices which 
behave like classical particles interacting via the (two-dimensional) Coulomb potential and submit¬ 
ted to a harmonic confining potential | 2121 1134112221122111213j . For extremely intense magnetic 
fields, the number of vortices tends to infinity, and the limit problem becomes that of Wigner’s 
crystallization (see Section 110)1) . as shown in |214j . This explains, although it has not been proved 
rigorously yet, why the hexagonal lattice appears in superconductors. We refer to [223 ] for a more 
detailed presentation of this problem. 

2.7.3. Ohta-Kawasaki model. The Ohta-Kawasaki model describes phase separation in copoly¬ 
mer systems [ 184l| . In its simplest version, it consists in minimizing the energy functional 

e|d{M= 1}| + i [ [ (u(x) - u)V(\x - y\){u{y) - u)dxdy 

1 J o J o 

where fl is a bounded domain of R d and where u is allowed to take only the two values ±1, each 
one corresponding to a phase of the system. The first term is the perimeter of the set where u = 1 
whereas the second term describes the interaction between the two phases. The potential V is often 
assumed to be the Coulomb or Yukawa interaction, with periodic boundary conditions on f 1 if it 
happens to be a cube. The question is to determine the optimal configurations, as the parameters 
e and u vary. If u = — 1, then the energy simplifies into 

e\dE\ + 2 f j V(\x~y\)dxdy (45) 

J e J E 

with E = {u = 1}. A regularized version of the model consists in minimizing the energy functional 

f (e' 2 \X7u(x)\ 2 + F(u(x))\dx+ f f (u(x) — u)V{\x — y\)(u(y) — u) dx dy, (46) 

Jn ^ ' Jn Jn 

where F is a non-negative function having as unique minimum points u = ±1. In the limit e —> 0, 
this problem becomes equivalent to (l45l) . 

In dimension one, it has been proved that the minimizer is periodic jl8l , IT] 1199; ! 53 i ; 252l ’ 1~12 j 
if u = 0, for all £ > 0. Very few results exist in higher dimension pTO] [54] . In the limit where one 
phase is strongly favored (u ~ —1) and £ —> 0, it has been proved that, here again, the opposite 
phase u = 1 is a solution to Wigner crystallization problem [1T71, [ 118] . 

A proof of crystallization on the hexagonal lattice (in 2D) has been recently given in [23] for a 
different copolymer model. In this theory, the second term in (IT51) is replaced by the Wasserstein 
distance W to the Lebesgue measure. Hence, the energy is defined for point measures y having 
their support in fl C R 2 . It reads 

£( v )=£ vV(( 4) + W(ln,M), 

2 :Esupp(/x) 




















The crystallization conjecture: a review 


23 


and crystallization is proved for any £ > 0 if 12 has appropriate symmetries, and if e is sufficiently 
small (or equivalently if e > 0 is fixed and |I2| —> oo). In the limit £ —> 0 it had been proved previ¬ 
ously that the hexagonal lattice minimizes the Wasserstein distance to the Lebesgue measure |182 j. 


2.8. The macroscopic object and its microscopic structure. We mentioned above the ques¬ 
tion of proving the existence of a macroscopic measure M, obtained as the weak limit (fT4l) . This is 
related to the formation of a macroscopic object. Another question is to know what kind of object 
is formed, that is, to compute the measure M. 

This problem seems different from the local behavior of the particles. However, the hexagon in 
Figure [T| indicates that a link exists with the microscopic scale. Indeed, if exact crystallization is 
assumed, that is, if the particles are restricted to be on the vertices of a periodic lattice for all N, 
then it is possible to write a limit minimization problem for the surface energy, which coincides with 
the second-order term in the development of E(N). This term is of order N ( d_1 )/ d . It has been 
proposed by Wulff j249j , and proved rigorously for a hard sphere mode{f| in dimension two by Au 
Yeung, Friesecke and Schmidt in |171 [220 ]. This work is based on results by Radin et al [ 1941 1 li) . 
We refer for instance to [25] [28] [52] for similar results on the Ising model. 

These works assume that the particles form a subset of a given periodic lattice for all N, which 
is true only for very specific interaction potentials V. It would be interesting to generalize these 
results to more general cases. However, this problem is a priori a very difficult one, since a good 
knowledge of the leading order term of E(N) is needed to understand the next one. And this is 
exactly the crystallization conjecture. 


3. Extensions 

3.1. Positive temperature. Until now we have only considered the problem of minimizing the 
energy, that is, we have assumed the temperature to be 0. As a matter of fact, it seems intuitive 
that crystallization only occurs for small temperature [ 196l| . At positive temperature T > 0, the 
problem is more complicated to state. We discuss here a possible formulation and mention some 
important references. 

The point particles should now be replaced by a probability density on R dJV that only describes 
random positions of the N particles. This probability is found by minimizing the free energy, which 
is the sum of the energy and of an entropy term multiplied by the temperature T. The entropy 
favors extended systems and a high value of T will spread out the probability density. In the whole 
space the entropy would make all the particles fly apart, and one therefore needs to confine them. 
This may be done as in Section 12.41 by imposing that the system is in a bounded domain 12 and that 
the volume |I2| is proportional to the (average) number of particles N. Another possibility would 
be to work in the whole space R d , and add a confining external potential as was done in CUl) . but 
we will (almost) not discuss this here. 


4 That is, V = +oo sur [0,1 — ej, V = 0 on [1 + £,oo) and min(l/) = V(l) < 0, with small e. 
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Formulation in the canonical ensemble. The probability that we have to consider is the Gibbs 
probability measure^ 


Pq,n,t{x i, Xn ) 


£n( x 1’-’ x n) 

e r- 


e T ax i • • • axN 

J n N 


(47) 


This distribution concentrates on the minima of £n when T —> 0. This probability measure is 
obtained by minimizing the free energy 


Fq(N,T ) := min 

P' symmetric measure 
on Q N , with 

m fn N p ‘ — 1 


{ M J N Xn ) P '( x l’ — > dx\ • • • ctejv 

+ -^r / ijv) logP'fe,■ ■ ■ dijv 1, (48) 

A ! J 


in which the first term is the energy of the system, and the second one is the opposite of the entropy. 
The symmetry of P' accounts for the fact that the particles are identical and indistinguishable. 
The Boltzmann coefficient 1 /N\ appearing in front of all the integrals is here to count each possible 
configuration only once. In (l48l) we can replace P' by the probability measure P = P' /(IV!) at the 
expense of an additive term Tlog(N\). At T = 0 an optimal P is a delta measure at one minimum 
of £/y But when T > 0, the entropy term imposes that P be absolutely continuous with respect to 
the Lebesgue measure. 

The solution of the minimization problem (1481) is unique, given by P' = (N\)P^^,t in (1771) . and 
satisfies 

Fq{N,T) = -T log Z n (N, T), 

where 

Z n (N,T) = l- [ e -^Gi ^ n )/ t dxi ■ ■ ■ dx N 
A ■ Jn N 

is called the (canonical) partition function. As before we look at the coupled limit AC —> oo with 
|f2jv| —> oo and AT/|f2jv| —> p > 0, a fixed parameter. Using the stability of V, £n > ACeoo, we find 


Fn N (A, T) > -T log ~ Ac( eoo + T(log p- 1)) 


by Stirling’s formula. Similarly, by Jensen’s inequality, 


Fn N (N,T) < -Tlog 


-1 Oat | N f n N£N/T\ 

[~ivr e N J 


NT (log p — 1) 


Np 


U(|x|) dx. 


A similar estimate exists when V is not integrable at 0. We deduce that Fq n (AC, T) behaves linearly 
in N. After a little more work (209 ]. one can actually prove that the thermodynamic limit exists 


f(P,T) 


lim 

N—t oo 

\n N \—>oo 
N/\h N \-+p 


Fq n (AC, T) 
N 


5 It is also possible to consider the Hamiltonian 0 instead of Sjy, but the Gibbs measure e 'H-n/t can then be 
factorized and the variables pi do not play any role. The situation is different in the quantum case fSection l3.3D . 
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as in C21) . 

In order to formalize the crystallization problem at positive temperature, it is convenient to 
consider the empirical measures (also called fc-point correlation functions |88| Il73| ). which are 
similar to the measure pjy introduced in (11511 . To be more precise, we define the family of measures, 
obtained by integrating with respect to all variables except k of them: 

<k) AT! f 

/TiJv.rOEii-i^fc) = , AT ‘ \ i / Pn,N,T(xi,...,Xk,yk+i,-,y N )dy k +i---dy N . (49) 

(iv — fcj! jQN-k +1 

(k) 

A natural definition of crystallization is that the measures / 1 N T locally converge, after extraction 
of a subsequence Nj and applying an appropriate translation t 3 to some locally finite measures in 
the thermodynamic limit, 

Nj,T\ ']> ^p,T 

where these measures are invariant under the action of a (maximal) lattice group G: 

VgeG, yW.fa + g, ...,x k + g) = p ( p ^(x 1 , ...,x k ). (50) 

Here again the limit will depend on the subsequence and on the sequence of domains (Hat). Different 
weak limits might have different invariance groups. The group G in (1501) is assumed to leave invariant 
all the possible weak limits up to rotations (for chosen T and p and all Nj and f2jv 3 - for which the 
limit exists). 

In general, when p( 2 \xi, x^) — p^\xi)p^\x 2 ) does not tend to zero when \x\ — a^l —> oo, one 
says that there is long range order, meaning that two particles far away are correlated, although 
the system is not necessarily periodic. Many works have been devoted to the proof that systems 
exhibit long range order, without reaching exact periodicity. 

For physical systems one might expect that there will be a periodic (or at least long-range) order 
for all T < T c (p) for a small critical temperature T c (p), and that the system will be translation- 
invariant, with unique equilibrium states, for T > T' c (p ) > T c (p). It is usually believed that the 
function f(p, T ) is piecewise real-analytic in (T, p) and that the curves of non-analyticity correspond 
to phase transitions IM . 

In the zero-temperature case, we only considered the measure p^ 1 ', which appears in (1161) . The 
reason is that, in such a case, a solution to (|48ll is 

Pn,N{yi,—,yN) = -jyi <W (1) (2 /i)'' ‘ (Vn) 

' tree N 

where x \,..., Xn is a solution to problem © . Because of this specific form, local convergence for 
k = 1 to p implies that all the other empirical measures p^ k \x 1 ,... ,x k ) automatically converge to 
p{x 1 ) • • • p(xk)- Such a property cannot hold in general for T > 0. In case of crystallization at 
T = 0, it is usually expected that the fc-particle densities p p !p should converge in the limit T —> 0 
towards the uniform average of translations and rotations of the crystal: 

k 

l im n Mp^(zi,---:Zfc) = A [ I f[p(R(xj-T))dRdr := p ik) 0 (xi,...,x k ). 

T ^° M Jq Jso(d)j-J[ 

Then p ^ is translation-invariant (in particular, p p jj is constant), but p^ 0 is periodic in X\ — X 2 - 
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Grand canonical ensemble. The problem is often stated in the grand-canonical ensemble, for 
which the algebra is simpler. This corresponds to assuming that the number of particles N in fl 
is also a random variable. It is then customary to consider the Laplace transform of the measures 
N i—^ Zq(N,T), introducing a variable p = ( T/N)\ogz dual to N, called the chemical potential: 

Zn(z,T) := £ e» N ' T Z n {N,T) = 1 + £ ^ f e~ e -' T . 

N>0 N>1 1 ' 


The series on the right is conveniently written in terms of the fugacity z = e ll / T . The stability 
£n > eoo-/V makes the sum convergent, with Zq(z,T) < exp(|f2|e^ _e ° 0 ^ ,/T ). By varying z > 0 (or 
p), one can obtain some information on the canonical case. The problem is then the same as before, 
except that we cannot divide by N which is an unknown variable, so we instead divide by |fl| and 
define 


f(z,T) 


lim 

|Q|—>-oo 


-T log Z n (z,T) 

\n\ 


which now satisfies —T|n|e^ e °°)/ T < f(z,T ) < 0. In a similar manner, we define the fc-particle 
grand-canonical correlation function by 


P , q'!z.t( x i> ■■■i x k) 


= Z n (z,T)- 1 z k e~ 


. 


E- 

n> 1 


y k-\-n 


n\ 


£(xj,.. x fc + n ) 

e ^ d,Xk -\-1 • • * dxk+ n I (^1) 


and ask the same questions as before about its possible local weak limits. 


Known results. There are many results showing the absence of crystallization in several situa¬ 
tions, but almost none proving its existence. The first classical theorems are due to Ruelle and 
Penrose (20911205 .. 189} who studied the convergence of the series at small fugacity z (hence, liigh- 
temperature for fixed p). They showed that the grand-canonical free energy as well as the correlation 
functions are all convergent series in the parameter z = e M / T , with a radius of convergence at least 
equal to 

/ /■ \ -i 


Rmin{T) = e^' T - X 


e -V{x)/T 


dx I 


(52) 


\JR d 


Furthermore, the correlation functions must be translation-invariant when 2 and T are in this range. 
Since eoo < 0 with strict inequality if min V < 0, the radius R m i n (T ) shrinks exponentially fast 
when T —> 0. The main tool here is a system of equations for the correlation functions called the 
Kirkwood-Salsburg equations , that allows to derive uniform bounds of the form (aq, ...,Xk) < C k 
at small fugacity. The bound on the radius of convergence was recently improved in •19.1 . The 
series of f(z,T) in terms of 2 is called the Mayer expansion. It is possible to invert the value 
of the constant density p(z,T) = ft 1 ' 1 in terms of 2 and to express / as a function of p and T. 
Since p(z, T) = z + 0(z 2 ), this is a low-density regime. The corresponding series (called the Virial 
expansion) was studied in 11511 . 

The properties of the system in the thermodynamic limit depend in a crucial way on the space 
dimension d, a phenomenon that is not present at T = 0. For one-dimensional systems with a 
fastly-decaying potential, a result due to van Hove in ( 237] (see also [209 '. Thm 5.6.7] and (207 , 67]) 
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states that f(z,T) is actually real-analytic for all z > 0 with unique correlation functions, hence 
there is no phase transition. Given that crystallization has been proved for a large variety of 
models in dimension one at zero temperature, this shows that the phase transition occurs exactly 
at T c (p) = 0. We notice that phase transitions can occur for an integrable potential decaying as 
slowly as l/|cc| 2 [74l[97] . A similar result called the Mermin-Wagner theorem [17811176111771 [Ml 196] 
exists in two dimensions for very short-range potentials, except that only the translation-invariance 
is known, and there does not seem to be any information on the range of analyticity. 

Given that there is no crystalline order in one and two-dimensional systems at positive temper¬ 
ature, but that crystallization is expected at T = 0, a natural question is to ask what is happening 
when T —> 0. Recent works [591 11311 fl33 j investigate the exponentially small region of convergence 
of the Mayer series close to T = p = 0, where things can be proved in any dimension. Namely, 
they considered the limit T —> 0 and p —» 0 with the constraint that Tlogp —> v (recall that 
p ~ z = at small activity, so this is similar to fixing p). They showed that, in any dimension, 
the system behaves as a gaz of finite clusters of particles which essentially do not interact. The 
number of particles K in the clusters is determined by the parameter v and it can be as big as we 
want if all the E(K) have minimizers. Each cluster of K particles solves the problem E(K). For 
K oo the cluster converges to the zero-temperature crystal configuration. 

Even without more information about the limiting states, it is possible to show in most situations 
that the correlation functions satisfy a pointwise bound of the form 

dn^ T ( x u—> x k) < (C z , T ) k (53) 

and therefore have local weak-* limits, after extraction of subsequences. The positive temperature 
makes the weak convergence of the correlation functions much easier than at T = 0 (see Section l2J2ll . 
but the constant Cfyr always blows up when T —> 0. The first simple case is that of a non¬ 
negative potential V: using in (ISlT) that £n{x\,...,X k+ n ) > £n(%Ii+u ■•;%k+n) we immediately 
obtain the pointwise bound (l53l) with C Z) t = z. The most general case was covered by Ruelle 
in [208 ] who proved (1 53 [I for super-stable interactions V. Super-stability means that there exists 
constants A, B, r > 0 such that 

v (\xi~Xj\)> (An(k) 2 - Bra(fc)) (54) 

1 <i<j<N fce z d 

where n(k) is the number of particles in the cube kr + [—r/2,r/2) d . Any continuous non-negative 
interaction with R(0) >0 satisfies (1541) for r small enough. Hence any stable interaction can be 
made super-stable by slightly increasing its values in a neighborhood of 0. Ruelle’s constant C Z) t 
is, however, a complicated function of the parameters. It was shown in [ IDI] that one can take 
C z ,t = ze ™(°)/ T when ui > 0. 

The Jellium problem described in Section 12.61 can be defined at positive temperature and it 
is an exception for which crystallization occurs for any T and any density p in dimension d = 
1 [14511391 17]. This is not surprising, since screening effects should make the effective potential 
be integrable at infinity, but probably not decay faster than 1/r 2 . Numerical simulations indicate 
that, in dimension d > 2, there exists a critical temperature T c > 0 such that, if T > T c , then 
Jellium is not crystallized (42[ 112511191111031 [9[[64[ [8l 12281 l69ll43l [30] . The link between a trapped 
Coulomb gas in the mean-field limit and the positive-temperature Jellium problem was recently 
studied in [1501 1149) , similarly to what we discussed in (l43l) . 
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Link with random matrices. The N eigenvalues of an N x N matrix with random coefficients 
are, in some situation, distributed according to the Gibbs measure of a gas of particles in an external 
trapping potential. The effective interaction will usually be V(x) = — log|x|, and the dimension 
d = 1 (if the eigenvalues are real) or d = 2 (if they are complex). If the entries of the matrix are 
independent Gaussian variables, the statistical distribution of the eigenvalues Ai,..., A at is given by 
the Gibbs measure (1471) with the mean-field energy 

1 N 
£ n (Xi,...,X n ) =-— Y log \X k - Ai| + ^ 1 Ay| 2 . 

i<k<e<N j =l 

For hermitian matrices ( GUE , that is, Gaussian Unitary Ensemble ), the problem is set in Q = 
R, since the eigenvalues have no imaginary part. In such a case, the temperature is equal to 
T = 1/(2 AT). If one imposes that the matrices have real coefficients ( GOE, that is, Gaussian 
Orthogonal Ensemble ), the temperature is T = 1 /TV. When considering complex matrices without 
any symmetry assumption (Ginibre ensemble), we have the same formula, but the A j are now in 
fl = C = K 2 and the temperature is T = 1/(2TV). 

It is also possible to consider unitary or orthogonal matrices ( CUE for Circular Unitary En¬ 
semble, and COE for Circular Orthogonal Ensemble, respectively), using the uniform law on this 
compact subset of matrices. Then, the eigenvalues are distributed according to the Gibbs mea¬ 
sure 63 with 

£jv(Ai, ..., Ajv) = —— Y, l°g|Afe-A^|, 

V 1 <k<i<N 

this time restricted to the unit circle f l = S 1 . 

Studying the eigenvalues of random matrices and the link with Coulomb gas is a very active 
subject, which started with the seminal works of Wigner [24711248] and Dyson [7D., [71;i7Z } [77 }lI75| . 
The interest in the set of matrices we just mentioned is that they allow for explicit computation of 
empirical measures, hence a good knowledge of the statistics of theses eigenvalues. Since T behaves 
like 1 /TV, the first order corresponds to the zero-temperature setting. The average distribution of 
the eigenvalues is given by the measure <r solution to the minimization problem in (1411) . The next 
order is more complex and its link with the crystal problem is less clear m We refer for instance 
to [[174 .', 92’. fl3] for a detailed study of the subject. 

3.2. Several types of particles. In order to deal with long-range interactions (for instance 
Coulomb potential), as in Wigner problem presented in Section [77771 it is possible to add a back¬ 
ground homogeneous density making the system globally neutral. Another model, more important 
from a practical viewpoint, is the case of two (or more) different types of atoms or ions, with dif¬ 
ferent charges. One can think for instance of sodium chloride crystal, which is made of two face 
centered cubic lattices, one (Na + ions) shifted with respect to the other (Cl” ions). 

For the sake of simplicity, let us consider only two types of particles. The interaction between 
two identical particles is different from the interaction between two different ones. We are thus led 
to the energy 

£ n \, n -2 (*dj • * * i , yii *■■) y n 2 ) 

JVi n 2 

= Y 1/11 (l Xi “ x il) + Y ^22 ( 1 ^ - Vj\) + Y Y Vl2 (l x i ~ Vi\) 

1<z< < 7<AT2 i= 1 j= 1 
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where Xj and yi are the positions of the particles of each type. We study the limit N\ , N2 —» 00 , 
possibly imposing a link between Ni and N 2 , accounting for a charge difference between the two 
types of particles. Thinking of a 3D crystal composed of charges of opposite sign q\ et —q 2 , we 
assume that 


9? 


Vu{\x\) ~ 

| x | —^00 \X\ 


<& 


MM)~ fl’ 

| x | —Yoo \X\ 


Vl 2 (\x\) 


~ gl <?2 

|x|—>oo Ixl 


and we impose the neutrality condition q\N\ — < 72-^2 —> 0 in the limit. For such a classical model, 
the Coulomb interaction is not adapted, since the energy tends to —00 as two particles of opposite 
charge get closer to each other, and the model is unstable (note, however, that it is stable in the 
quantum case, as it was proved by Dyson-Lenard |75], 1 152] and Lieb-Thirring | 165j ). Hence, one 
needs to assume that the potentials Vn V 22 and V 12 are repulsive at short distance \x\ —» 0. The 
simplest choice is to take 


Vij(\x\) = 2(Sij - l/2)q i q j W(x), with W(x) = -—-. 

\x\ 

Since W has a positive Fourier transform, the same calculation as in m gives that the total 
interaction is bounded from below by 


^n±,n 2 (xi ? ■ • ■ ? xn-^ , y \, ■ ■., yN 2 ) 

1 

> 


2(2n) d / 2 


r _ 

JVi n 2 

/ W(k) 

< 7 l ^ e ik ' x i - q 2 Yl e ik ' Vk 

m d 

i= 1 fe= 1 


N 1 qf + N 2 ql N 1 q\ + N 2 ql 
ak --- y, > --- y, 


and the model is stable. 

Several conjectures have been made concerning the optimal lattices m, but we do not know 
any result on the crystal problem with several types of particles. Thinking of crystalline structures 
currently observed in nature, it is a highly important question from a physical viewpoint. The 
existence of the thermodynamic limit and estimates on the correlation functions were proved in [ 101 J. 
A review of known results in 3D for high temperature (hence without crystallization) is given in 
[43] . In [195 ] Radin considers special short-range potentials for two types of particles, and proves 
that crystallization fails, but the minimizers are quasi-periodic. 


3.3. Quantum models. In the classical models studied so far, the kinetic energy of the particles 
does not play any role, since we deal with minimizers or Gibbs states. The term Y^jLi |fb| 2 /(2m) 
in © disappears in the minimization problem, and factors out and gives a Gaussian at positive 
temperature in dUD. The situation is different in quantum mechanics, in which there is a link 
between velocity and position, in order to respect Heisenberg’s uncertainty principle. This makes 
the kinetic energy dependent on the positions of the particles. More precisely, quantum mechan¬ 
ics principles imply that pj should be replaced by the differential operator —ihS7 Xj and that the 
Hamiltonian ...,pat,Xi, ... ,Xn ) in (HD should be replaced by the differential operator 

N n 2 

h n = ~y~]+ X] v ( Xk - x *)• ( 55 ) 

j=1 1 <k<e<N 
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This operator acts on L 2 (Q N ), where = R d for an unconfined system, and where SI is a bounded 
domain if the system is confined (with suitable boundary conditions). Since the particles are 
indistinguishable, we work with a subspace of L 2 (il N ) consisting of functions having a prescribed 
symmetry property. In nature one can find two types of particles: bosons and fermions. For 
bosons, we use the subspace L 2 (fl N ) of functions which are symmetric with respect to variable 
permutations. For fermions, we use the subspace L 2 (Cl N ) of functions which are antisymmetric. 
Properties of the system in the limit TV —> oo depend on the chosen symmetry class. For the sake 
of simplicity, we asume that the particles have no spin. 

The classical problems studied so far read, in the quantum case, 


E a / s ,n(N) 


inf „ [ i^—\V^{xi,...,x N )\ 2 + V V(x k -xe)\'^(x 1 ,...,x N )\' 2 \dxi---dx N , 

*eL2/. ( nVn" l 2m ^ \ 


/i*r=i 

for the energy at T = 0, 


(56) 


F a /s,n,T(N) = -T log (tr^^e 
for the canonical free energy, and 


F a /s,n,z,T{N) — —T log 


,N>0 


^tr L 2 (oJV) e Hn/t 


for the grand-canonical free energy. Note the absence of the Boltzmann factor, due to the restriction 
that we work in the subspaces L 2 ^ s {Cl N ). If we instead work in the full space L 2 (fl N ) (indistin¬ 
guishable particles sometimes called “boltzons”), then we need to put back the 1/TV! coefficient. In 
the semi-classical limit h —> 0, these (free) energies converge (up to a constant which diverges like 
log H at positive temperature) to the corresponding classical (free) energies. 

Quantum mechanics is by nature a probabilistic theory and the study of crystallization relies on 
the weak limit of the empirical measures, as in the case of the positive temperature classical model. 
For instance, at zero temperature one can study the limit of the fc-particle densities 




TV! 


(TV - k)\ J n s-H 


I'Ffi.ArOci, ...,x k ,y k +i,—,yN )\ 2 dy k + i ■ ■ • dy N 


where is a minimizer of problem (1561) (this minimizer is always unique, up to a phase, for 

bosons, but it is not necessarily unique for fermions). A similar formula is valid at positive tern- 
perature but we do not provide the details. The measure /.4 q n corresponds to the one defined in 
classical mechanics and it does not carry all the information on the quantum system. It is more 
relevant to study the limit of k-particle density operators, which are defined by their integral kernel 

(k) / / / \ 

r iVl,N\ X ^i "•> X ki x li X k ) 

TV! r 


(N-k)lJ nN - k 


'i’n,N(x 1 ,...,x k ,y k+1 ,...,y N )'!>n,N{x , 1 , ...,x’ k ,y k+ 1, ...,2/ at ) dy k +i ■ ■■dy N , (57) 
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(k) 

and where the diagonal part coincides with N . We say that the system crystallizes if theses 
operators locally converge to operators 7 ^) which commute with translations of a (maximal) lattice 
G, that is, 

Vg e G, 7^(27 +g , ...,x k +g,x [ + g,...,x' k + g) = 'y^ N (x 1 , ...,x k ,x[, ...,x k ). 

If V is a stable potential, then the quantum mechanical problem is also stable since the Laplacian 
is a non-negative operator, and the thermodynamic limit can be shown to exist mg. However, 
potentials V that are unstable classically can become stable in quantum mechanics. This the case 
for the Coulomb potential in dimension d = 3 with two kinds of particles of opposite charge which 
has been shown to be stable for fermions [1651116611161] but not for bosons [731 ISUl [1631116411227) . 
The existence of the thermodynamic limit in the Coulomb case is proved in various situations |1561 
[941 l84l 1121111221 i27| . Correlation functions were studied in this case in (101 , 94] , 

Few results have been proved for the crystal problem in the case of continuum quantum systems. 
In particular, one could think that, when a classical system exhibits crystallization, so does the 
quantum corresponding system if h is sufficiently small. This has not been studied, to our knowledge, 
except in the case of Coulomb gas (quantum Jellium) for which Kunz [145] and Brascamp-Lieb [32] 
have proved crystallization for small density p in dimension d = 1. After a change of scale, assuming 
that p is small is equivalent to assuming K is small, so the situation is indeed a semi-classical limit. 
Crystallization for ID quantum Jellium (at any density p and any temperature T) has been recently 
proved by Jansen and Jung 132| . 

Many of the results on the non-existence of crystals in the classical case have been extended 
to quantum systems. Using the Feynmann-Kac formula (expressing all the quantum objects in 
terms of paths in the classical problem), Ginibre proved in (107 . 10811109] that the grand-canonical 
free energy as well as the fc-particle density matrices are convergent series in terms of the fugacity 
z. Hence there is no phase transition for z small enough and the fc-particle density matrices are 
translation-invariant. The uniform bounds (15311 allowing to define local weak limits have been 
generalized to quantum systems in ( 1011 (94l 182 , 1188J . 

Some results have been proved for quantum systems described by nonlinear models, such as 
Thomas-Fermi or Hartree. In the Coulomb case, assuming that the nuclei are classical particles 
with positive charge and are distributed on a lattice, it has been proved for convex models that 
the electrons are periodically arranged [1621 1551 [501 ISI1I37|. If in addition one optimizes over the 
positions of the nuclei, then crystallization is only known in ID for Thomas-Fermi type models |26| . 

For a general interaction, one could study the Hartree problem that is the mean-field approxi¬ 
mation of the bosonic many-particle system (153) . It relies on the nonlinear energy functional 

£ ( u )= 7T- I \Vu(x)\ 2 dx+ / f V{\x-y\)\u{x)\ 2 \u(y)\ 2 dxdy, 

1/ ^ 1/ v 

with f \u(x)\ 2 dx = N, (58) 

J n 

which is similar in spirit to (1461) . The numerical simulation presented in Figure [Til shows the 
occurrence of crystallization for this model in ID for a Lennard-Jones interaction. Similar results 
have been observed for the potential V(r) = 1 ( 0 , 1 ) C 7 ') i n higher dimensions [192 , 13511136 ; :4lf51 2451 

n~ 44 i n~ 3 rn n~ 2 ] . 
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3.4. Discrete systems. In our review, we focused on continuous systems, defined in the whole 
space or in a domain fljv which grows as N —> +oo. An important literature is devoted to the 
study of discrete systems. Such systems are defined on a lattice G C R d , without assuming a priori 
that the corresponding states are G-periodic. We expect that proving the occurrence of periodicity 
is slightly easier, since a natural periodic lattice is already present in the definition of the system. 
Many rigorous results have been proved for this kind of problems, although important questions 
are still unsolved. The models considered can be either quantum or classical. 

Examples of such systems are the (classical or quantum) Heisenberg or Ising models. Two main 
regimes are usually dealt with: the ferromagnetic one, in which spins tend to be aligned with each 
other, and antiferromagnetic in which spins are preferably of alternate sign. In this latter case, 
crystallization gives a periodic lattice which size is twice that of the original one. 

In 1986, Kennedy and Lieb have considered two systems of this type. In [1371 113811155] . they 
study electrons on a lattice, submitted to a pointwise interaction with fixed particles of opposite 
spin. They prove that the electrons are located on a sub-lattice. In jI3U] they consider a ID system 
on the lattice Z. This model describes for instance deformations of a polyacetylene molecule. They 
prove that the minimizer is periodic of period 2, a phenomenon called Peierls instability. This 
result has been further developed in |158[ |157[ |159j. It has been generalized in m, and extended 
to the hexagonal lattice in 2D in [93] . 

Apart from systems with analytical solutions, an important method for studying classical or 
quantum spin systems is the reflection positivity method. This strategy has been introduced in 
field theory [187] . then adapted and developed in the case of spin systems [102 , 76, lOOL [68 : i66]. 
This method aims at proving phase transitions and long-range order. However, it does not always 
allow to conclude that the system is periodic. For recent examples of application of this theory to 
crystallization problems, see for instance mmmm- 


Conclusion 

We have described several aspects of an important problem arising in physics and which, in spite 
of an intense activity, is still not completely understood mathematically. In addition to the famous 
crystallization problem, several questions have been mentioned, some of which are probably more 
at hand than others. Some progress in any of these directions would be of high interest and would 
improve the theoretical understanding of the structure of matter at the microscopic scale. We hope 
that this article will stimulate further research in these directions. 
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Figure 5. Left: &(S, s) as a function of s for the square lattice (S = I 2 ) and the hexagonal one (S = Shex)- 
Right: the relative difference (£ 2 ( 12 , s) — C 2 (Shex, s))/|C 2 (Shex; s)|. It shows that the hexagonal lattice 
energy is lower than that of the square lattice for all s > 0, as it is proved in 1981 l48l 17911651 . 




Figure 6. Left: 62 (a) as a function of a for the square lattice (S = I 2 ) and the hexagonal one (S = She x)- 
Right: the relative difference ( 62 ( 12 , a) — ^(Shex, a:))/6*2(5'hex, a). It shows that the hexagonal lattice energy 
is lower than that of the square lattice for all s > 0, as it is proved in |l80 i. 
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Figure 7. Left: <^3 (S', s) as a function of s for different lattices ; FCC and BCC have energies which are very 
close to each other. Right: the relative difference (( 3 (BCC, s) — ( 3 (FCC, s)J/ |( 3 (FCC, s)|. It indicates that 
BCC should be the minimizer for 0 < s < 3/2, while FCC should be for s > 3/2. Proving this is still an 
open problem. The relative difference is of order 10~ 4 . 




Figure 8 . Left: # 3 ( 0 ) as a function of a for several lattices ; FCC and BCC have energies which are very 
close to each other. Right: the relative difference (# 3 (BCC, a) — 03 (FCC, a))/^(FCC, a). It indicates that 
FCC is the minimizer for a > 1, while BCC is for a < 1. Proving this is an open problem, which does not 
seem to imply the above result on the function ( 3 . 
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Figure 9. A colloidosome is a spherical shape made of colloids (polystyren molecules here), which is described 
by Thomson’s model. This model consists in minimizing the interaction |x| _s for particles on the sphere 
S 2 Cl 3 . Scanning microscope picture, from (66| . © AAAS. 



Figure 10. Left: Experimental pictures of fast rotating Bose-Einstein condensates: the number of vortices 
increases with the rotation velocity. The experiments have been conducted by Ketterle’s team |Tj at 
MIT in 2001. © AAAS. Right: Numerical simulation of the Gross-Pitaevskii equation with the software 
GPELab [mill, reproducing the vortex lattice in the corresponding regime. 
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Figure 11. Numerical calculation of the minimizing density |u| 2 of the Hartree energy (1581) in ID for the 
truncated Lennard-Jones potential w(x) = min(10 3 , |ie| -12 —1*| _6 ). Here N = 10, h 2 /(2m) = 1, the average 
density is p = 1 and we have used Dirichlet boundary conditions, on a mesh with 2 x 10 3 points on the 
interval [0,10]. 



